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ABSTRACT 


Geophysical  models  constitute  an  important  component  of  calibration  for  nuclear  explosion  monitoring.  In 
order  to  keep  them  relevant  to  current  monitoring  problems,  we  have  enhanced  the  models  in  several 
different  ways.  We  will  focus  on  discussing  several  of  the  enhancements  here.  We  have  significantly 
improved  upon  our  surface  wave  model  by  expanding  the  region  of  the  existing  model  of  Eurasia  and  North 
Africa,  south  to  cover  all  of  Africa  and  north  to  cover  the  polar  region  into  Alaska  and  Canada.  We  have 
also  improved  on  our  coverage  in  existing  regions  by  including  station- station  dispersion  paths  based  on 
ambient  seismic  noise.  Surface  waves  are  also  being  used  to  derive  the  crust  and  upper  mantle  velocity 
structure  of  these  regions,  including  important  parameters  such  as  crustal  thickness,  upper  mantle  velocity, 
and  lithospheric  thickness.  We  will  be  comparing  some  of  these  results  to  other  estimates  of  these 
parameters.  We  demonstrate  some  of  the  latest  enhancements  that  we  have  made  regarding  stochastic 
models.  For  instance,  we  have  been  using  the  Markov  Chain  Monte  Carlo  (MCMC)  technique  to  produce 
stochastic  models  in  the  Yellow  Sea  -  Korean  Peninsula  (YSKP)  region.  By  including  more  data  sets  (Love 
wave  dispersion  curves,  gravity)  and  more  data  from  the  existing  data  sets  (more  travel  times,  more 
receiver  functions,  more  Rayleigh  wave  dispersion  curves),  we  have  improved  the  lateral  resolution  of  the 
model  from  2  degrees  to  1  degree.  We  will  discuss  some  details  of  the  methodology  as  well  as  features  of 
the  model.  Finally,  we  will  focus  the  last  section  on  research  to  move  beyond  past  and  current  1-D,  2-D, 
and  2/4-D  methods  and  discuss  some  of  the  ongoing  efforts  to  transition  to  fully  3-D  models. 
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OBJECTIVE 


Geophysical  models  constitute  an  important  component  of  calibration  for  nuclear  explosion  monitoring. 

The  objective  of  any  geophysical  model  is  to  improve  prediction  of  travel  times,  amplitudes,  and  other 
observables  of  local,  regional,  and  teleseismic  phases,  in  order  to  improve  location,  detection,  association, 
identification  and  (in  the  event  of  an  explosion)  yield  determination.  In  this  paper,  we  will  focus  on  recent 
enhancements  to  many  of  the  models,  which  have  kept  them  up  to  date  and  relevant  for  current  monitoring 
problems.  We  first  discuss  improvements  to  the  surface  wave  models  which  have  expanded  to  cover  over 
half  of  the  globe  and  have  increased  in  resolution  due  to  the  increase  in  measurements  and  available 
seismic  data,  as  well  as  techniques  such  as  ambient  noise  correlation  that  allow  us  to  sample  new  paths. 
Tomographic  images  over  a  wide  period  range  are  robust  and  reliable  enough  to  allow  the  inversion  of 
highly  developed  models  that  cover  both  shallow  sediment  structure  and  deep  lithospheric  models,  even  in 
aseismic  regions.  The  results  can  also  be  focused  on  regional  or  station- specific  applications.  Geophysical 
models  such  as  the  Western  Eurasia  and  North  Africa  (WEN A)  model  can  be  updated  to  reflect  our 
increased  understanding  of  regions.  They  continue  to  both  be  useful  in  their  own  right,  and  also  serve  as 
critical  starting  models  for  more  data  driven  models,  such  as  the  Seismic  Location  Baseline  Model  (SLBM) 
(Myers  et  al.,  these  Proceedings).  Finally,  stochastic  models,  which  require  the  most  sophisticated  inversion 
technique  and  computational  power,  are  increasingly  becoming  easier  to  implement,  allowing  broader 
applications  of  the  method.  Future  enhancements  to  geophysical  models  will  increasingly  focus  on 
improving  resolution  and  combining  individual  data  sets  to  produce  a  fully  3-D  model  which  is  predictive 
of  a  range  of  observables. 

RESEARCH  ACCOMPLISHED 


Surface  Wave  Models 

Over  the  past  several  years,  Lawrence  Livermore  National  Laboratory  (LLNL)  has  been  developing  surface 
wave  models  for  nuclear  explosion  monitoring  (Pasyanos  et  al.,  2001;  Pasyanos,  2005).  Dispersion 
measurements  are  made  using  multiple  narrow-band  filters  on  deconvolved  displacement  data  from  the 
LLNL  Seismic  Research  Knowledge  Base  (SRKB).  We  continue  to  improve  upon  our  surface  wave  model 
by  adding  more  paths  by  taking  advantage  of  new  data  sets  and  also  by  revisiting  stations  with  more  recent 
events.  Most  recently,  we  have  added  measurements  from  stations  in  Alaska  and  Canada,  Siberia,  southern 
and  central  Africa,  and  southeast  Asia.  To  date,  over  120,000  seismograms  have  been  analyzed  to 
determine  the  individual  group  velocities  of  7-150  second  Rayleigh  and  Love  waves.  Overall,  we  have 
made  good  quality  dispersion  measurements  for  41,000  Rayleigh  and  30,000  Love  wave  paths.  Using  a 
conjugate  gradient  method,  we  then  tomographically  invert  these  measurements  to  produce  group  velocity 
maps  for  Love  and  Rayleigh  waves. 

The  group  velocity  models  have  been  enhanced  in  several  ways.  First,  with  more  measurements,  we  have 
been  able  to  expand  the  region  of  coverage  from  Eurasia  and  North  Africa  both  south  to  include  all  of 
Africa  and  north  (over  the  pole)  to  include  the  northern  parts  of  North  America  (Figure  1  and  Figure  2).  In 
these  extended  regions,  path  coverage  (in  particular  crossing  ray  paths)  is  sufficient  for  intermediate 
periods  where  we  have  many  measurements,  but  we  currently  do  not  have  sufficient  coverage  for  high 
resolution  at  shorter  periods  where  some  features  are  smeared.  We  anticipate  improvement  at  these  periods 
with  additional  dispersion  measurements.  Path  coverage  throughout  the  broad  region  has  been  further 
improved  by  the  use  of  cross-correlation  of  ambient  seismic  noise  to  derive  the  Green  function  between 
two  stations,  and  from  which  the  dispersion  characteristics  can  be  derived  (Ritzwoller  et  al.,  these 
Proceedings).  The  benefits  of  this  method  in  seismology  have  been  dramatically  demonstrated  for  southern 
California  in  Shapiro  et  al.  (2005).  Path  maps  for  Europe  (from  Yang  et  al.,  2007)  are  illustrated  in 
Figure  3. 
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20.0  second  Rayleigh 


20.0  second  Rayleigh 
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Figure  1.  Surface  wave  dispersion  maps  for  20  second  Rayleigh  waves,  (a)  Path  map  showing 

coverage  of  the  Eurasian  and  African  continents,  (b)  Group  velocity  map  of  the  same 
region  (Pasyanos,  2005). 


80.0  second  Rayleigh 
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Figure  2.  Surface  wave  dispersion  maps  for  80  second  Rayleigh  waves,  (a)  Path  map  showing 

coverage  of  the  Eurasian  and  African  continents,  (b)  Group  velocity  map  of  the  same 
region  (Pasyanos,  2005). 
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Noise-based  Dispersion  Paths 


Figure  3.  Cross-correlation  path  map  for  Europe.  Stations  are  indicated  by  green  triangles 
(measurements  from  Yang  et  al.,  2007). 

Velocity  Models 

One  application  of  the  surface  wave  models  is  to  invert  the  dispersion  results  to  derive  models  of  the  crust 
and  upper  mantle  structure.  This  is  particularly  useful  in  aseismic  regions  that  are  poorly  sampled  by  other 
data  sets.  Recently,  we  have  used  surface  wave  dispersion  maps  to  construct  a  3-D  lithospheric  velocity 
model  of  Africa  and  Arabia  (Pasyanos  and  Nyblade,  2007),  most  of  which  is  aseismic  and  not 
well-sampled  by  body  waves.  A  crustal  thickness  map  of  the  region  is  shown  in  Figure  4a.  In  North  Africa, 
there  is  a  crustal  thickness  contrast  between  the  West  Africa  Craton  and  the  East  Saharan  Shield.  There  is 
crustal  thinning  along  the  African  Rift  Systems  (West  and  Central,  as  well  as  the  East  African  Rift).  In 
central  and  southern  Africa,  crust  is  consistently  thick,  although  it  appears  to  be  thickest  under  the  Congo 
Craton.  Figure  4b  shows  a  slice  through  the  model  at  200  km  showing  the  faster  velocities  under  the  West 
African,  Congo,  and  Kalahari  Cratons  in  northern,  central,  and  southern  Africa,  respectively. 
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Crustal  thickness 


200  km  slice 


Figure  4.  a)  Crustal  thickness  map  of  Africa,  b)  S-wave  velocity  variation  in  the  mantle  under  Africa 
at  200  km  depth  (from  Pasyanos  and  Nyblade,  2007). 

By  applying  the  inversion  of  surface  waves  for  structure  to  broad  regions  we  can  start  to  map  out  important 
tectonic  features.  For  example,  the  velocities  of  long  period  surface  waves  are  sensitive  to  the  thickness  of 
the  lithosphere,  which  have  been  inverted  for  in  Figure  5.  Slower  group  velocities  at  short  periods  indicate 
the  lack  of  a  thick  lithospheric  lid  and  the  presence  of  a  slower,  asthenospheric  mantle  not  far  below  the 
crust.  These  are  found  at  plate  boundaries  whether  subducting  convergent  (Mediterranean,  Japan),  orogenic 
convergent  (Turkish-Iranian  Plateau,  Himalayas),  or  divergent  (Red  Sea,  oceanic  rifts).  Fast  velocities  at 
long  periods  indicate  the  presence  of  a  thick  lithospheric  lid,  as  well  illustrated  by  the  lithospheric  keels 
found  under  platforms  and  shields  such  as  the  West  African  Craton,  Congo  Craton,  Indian  Shield,  South 
China  Craton,  Baltic  Shield,  Ukrainian  Shield,  and  Russian  Platform. 
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Lithospheric  thickness 
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Figure  5.  Lithospheric  thickness  map  of  Eurasia  and  North  Africa.  Plate  boundaries  are  indicated  by 
the  solid  black  lines. 

Regional  Applications 

We  have  used  the  surface  wave  models,  mostly  in  conjunction  with  other  data  sets  such  as  receiver 
functions,  to  produce  regional  velocity  models.  Recent  studies  have  covered  Arabia  (Tkalcic  et  ah,  2006), 
Kenya  (Benoit  et  ah,  2006),  Kuwait  (Pasyanos  et  ah,  2007),  eastern  Turkey  (Gok  et  ah,  2007a),  eastern 
Mediterranean  (DiLuccio  and  Pasyanos,  2007),  Iraq  (Gok  et  ah,  2007b),  India  and  South  Africa. 

As  an  example,  the  KUW1  model  was  derived  from  the  joint  inversion  of  receiver  functions  and  surface 
wave  group  velocities  at  station  KUW  in  Kuwait.  It  is  a  45  km  thick  crust  with  8  km  of  sediments  and  very 
slow  upper  mantle  velocities.  Figure  6  shows  how  the  KUW1  model  (indicated  by  the  solid  blue  line) 
compares  to  other  models.  In  particular,  we  compare  this  velocity  model  to  models  determined  from  the 
simultaneous  relocation  of  local  seismic  events  using  VELEST.  These  VELEST  models  do  not  have 
sensitivity  in  the  shallow  sediments  and  deep  crust,  but  the  resulting  velocity  profiles  (shown  by  the  red, 
brown,  and  black  colored  dashed  lines)  are  very  similar  to  the  KUW  1  model  in  the  middle  crust,  regardless 
of  various  starting  models,  indicating  that  the  KUW1  is  also  applicable  for  locating  regional  events.  The 
preferred  VELEST  model  (where  the  sedimentary  column  and  deep  crust  are  set  to  KUW1)  is  indicated  by 
the  solid  purple  line.  The  resulting  models  have  significantly  thicker  crust  than  the  current  model  used  by 
the  Kuwait  National  Seismic  Network  (KNSN),  shown  in  green. 
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P-  wave  Model 


Velocity  (km/s) 

Figure  6.  P-wave  velocity  profiles  for  Kuwait.  See  text  for  details  (Pasyanos  et  al.,  2007). 

While  the  model  can  simultaneously  fit  surface  waves,  receiver  functions,  and  local  travel  times,  the  model 
does  not  appear  appropriate  for  regions  not  too  distant  from  the  station.  Figure  7  shows  the  travel  times  of 
local  Kuwaiti  events  and  regional  phases  (Pn,  Pg,  Sn,  Lg)  from  events  in  Iran  recorded  at  station  KUW, 
paths  which  primarily  sample  the  Zagros  Mts.  The  model  fits  reasonably  well  the  local  events 
(<1°  distance)  and  the  near  regional  events  (3-4°),  but  the  far  regional  events  which  have  longer  paths 
through  the  Zagros  require  faster  upper  mantle  P-wave  velocities  and  slower  crustal  P-wave  and  S-wave 
velocities. 


Figure  7.  Travel  times  from  local  events  and  regional  events  in  the  Zagros  recorded  at  station  KUW 
in  Kuwait  (Pasyanos  et  al.,  2007). 
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A  Priori  Models 

The  WENA  model  is  an  a  priori  three-dimensional  model  for  Western  Eurasia  and  North  Africa 
(Pasyanos  et  ah,  2004).  Models  like  WENA  can  serve  as  background  values  for  travel  time  correction 
surfaces  and  other  derived  parameters.  In  the  original  publication,  we  demonstrated  our  ability  of  the 
WENA  model  to  predict  a  variety  of  data  sets,  including  surface  wave  dispersion  maps  and  matching 
crustal  thickness  estimates  from  receiver  functions.  In  a  followup  paper,  it  was  also  demonstrated  that  it  is 
also  possible  to  improve  regional  travel-time  prediction  and  seismic  event  location  accuracy  using  this 
model  (Flanagan  et  al.,  2006).  The  authors  conclude  that  a  priori  models  are  directly  applicable  where  data 
coverage  limits  tomographic  and  empirical  approaches,  and  the  development  of  the  uncertainty  model 
enables  merging  of  a  priori  and  data-driven  approaches  using  Bayesian  techniques. 

This  model  continues  to  evolve  in  several  ways.  First,  the  WENA  model  was  joined  with  a  Los  Alamos 
model  for  Eastern  Eurasia  to  produce  a  single  model  for  Eurasia  and  North  Africa  (dubbed  the  Unified 
model).  Secondly,  crustal  thicknesses  in  the  model  have  been  replaced  by  a  1 -degree  crustal  thickness 
model  derived  from  kriging  thousands  of  individual  crustal  thickness  estimates.  The  model  continues  to  be 
tested  and  validated,  recently  to  provide  the  best  1-D  velocity  profiles  for  waveform  modeling. 

Most  recently,  the  attenuation  structure  of  the  model  is  being  updated  to  replace  the  original  placeholder 
values  that  were  provided  for  crystalline  crust  and  uppermost  mantle  layers.  These  updates  will  better 
reflect  the  tectonic  history  of  the  geophysical  province  (e.g.,  orogenic  zone  vs.  shield)  and  allow  the 
attenuation  predictions  of  the  model  to  more  closely  match  those  from  empirical  attenuation  studies,  such 
as  those  derived  from  Lg  coda  (e.g.,  Mitchell  et  al.,  1997). 

Stochastic  Geophysical  Models 

In  recent  years,  LLNL  has  been  developing  and  refining  the  use  of  the  Markov  Chain  Monte  Carlo 
(MCMC)  method  to  generate  three-dimensional,  data-driven  stochastic  earth  models.  We  use  MCMC  to 
sample  models  from  a  prior  distribution,  test  them  against  multiple  data  types  in  a  staged  approach,  and 
develop  a  posterior  distribution  of  seismic  models  that  are  most  consistent  with  multiple  seismic  data  sets 
(Pasyanos  et  al.,  2006).  This  approach  has  several  advantages  over  a  single  deterministic  model.  First,  we 
are  able  to  easily  incorporate  prior  information  on  the  model,  such  as  the  a  priori  geophysical  models  that 
were  discussed  earlier.  Secondly,  with  this  technique,  we  are  able  to  reconcile  different  data  types 
(such  as  body  waves  and  surface  waves)  that  can  be  used  to  constrain  the  model.  We  can  also  estimate  the 
uncertainties  of  model  parameters,  properly  migrating  data  uncertainties  into  model  uncertainties.  Finally, 
we  can  estimate  uncertainties  on  predicted  observable  signals,  such  as  would  be  required  to  apply  this 
model  as  a  correction  surface. 

While  this  approach  has  several  advantages  over  deterministic  models,  the  application  of  the  method  can  be 
computationally  expensive,  allowing  model  development  only  in  high-performance  computing 
environments.  For  example,  the  application  of  this  method  to  determine  the  crust  and  upper  mantle 
structure  of  the  Yellow  Sea  and  Korean  Peninsula  (YSKP)  region  using  surface  wave  dispersion  data,  body 
wave  travel  time  data,  gravity,  and  receiver  functions  has  necessitated  the  use  of  Livermore  Computing 
(LC)  machines  like  Thunder,  a  23  tflops  peak-performance  1024  Node  cluster  where  all  nodes  are  4-CPU, 
1.4  GHz  Itanium  2  "Madison  Tiger4"  nodes  (see  Figure  8).  As  we  wish  to  apply  this  methodology  to  even 
broader  regions  with  significantly  larger  data  sets,  one  key  objective  has  been  to  increase  the  efficiency  of 
our  runs.  Over  the  past  six  months,  we  have  solved  some  of  the  technical  problems  with  increasingly  larger 
data  sets  and  have  sped  up  the  calculations  significantly,  allowing  us  to  apply  this  to  bigger  problems. 
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Figure  8.  Livermore  Computing  (LC)  cluster  Thunder. 

Stochastic  models  represent  one  way  of  transitioning  to  3D  models,  although  there  are  a  number  of 
different  approaches.  Methods  that  developed  separately  are  starting  to  converge.  Surface  wave  models  are 
increasing  in  resolution,  allowing  detailed  shallow  structure  models.  Body  wave  models  have  increasingly 
sophisticated  modeling  of  the  ray  path,  allowing  more  data  to  constrain  the  model.  Data  sets  and 
computational  power  are  both  growing  rapidly.  On  June  6-7,  2007,  in  Berkeley,  CA,  LLNL  hosted  the 
“Workshop  on  Multi-Resolution  3D  Earth  Models  to  Predict  Key  Observables  in  Seismic  Monitoring  and 
Related  Fields”  (Harris  et  al.,  these  Proceedings).  Questions  raised  during  the  workshop,  such  as  model 
representation,  data  used  to  construct  the  model,  operational  constraints,  uncertainties,  etc.,  all  have  a 
significant  effect  on  the  best  approach. 

CONCLUSIONS  AND  RECOMMENDATIONS 


Geophysical  models  play  an  important  role  in  nuclear  explosion  monitoring,  allowing  us  to  capture  the 
effect  of  the  earth  on  seismic  signals  in  a  compact  and  self-consistent  manner.  As  the  availability  of  seismic 
data  continues  to  expand  and  computational  power  rapidly  grows,  we  have  increasing  power  to  truly  image 
the  high-resolution  features  of  the  earth  that  are  necessary  to  support  nuclear  explosion  monitoring 
operations.  We  have  presented  a  variety  of  models  currently  being  developed  at  LLNL.  As  data  availability 
and  computational  power  increase  the  models  can  be  enhanced  and  combined  in  an  effort  to  produce  fully 
3-D  models. 
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